Technical  Document  1226 

January  1988 

Daylight  and  Twilight  Sky 
Radiance  and  Terrestrial 

Irradiance 

C.H.  Waldman 


Titan  Systems,  Inc. 


The  views  and  conclusions  contained  in  this 
report  are  those  of  the  authors  and  should  not 
be  interpreted  as  representing  the  official 
policies  either  expressed  or  implied,  of  the 
Naval  Ocean  Systems  Center  or  the  U  S 
Government 


NOSCTD  1226 


NAVAL  OCEAN  SYSTEMS  CENTER 

San  Diego,  California  92152-5000 

R.  M.  HILLYER 
TeCiii  tical  Director 


E.  G.  SCHWEIZER,  CAPT.  USN 
Commander 


ADMINISTRATIVE  INFORMATION 

This  report  was  prepared  by  Titan  Systems.  Inc.  under  contract  N6600 1-84-D-0009  for  the  Naval 
Ocean  Systems  Center,  Code  843. 


Released  by 
J.P.  Meenan,  Head 
Electro-Optic  Devices 
Branch 


Under  authority  of 
M.S.  Kvigne,  Head 
Space  Systems  Division 


JJ 


UNCLASSIFIED 

SECURITY  CLASSIFICATION  OF  THIS  PAGE 


la  REPORT  SECURITY  CLASSIFICATION 

UNCLASSIFIED _ 

in.’securiiTY  classification  authority 


REPORT  DOCUMENTATION  PAGE 


lb.  RESTRICTIVE  MARKINQ8 


3  DISTRIBUTION/AVAll.AWUTY  OF  REPORT 


2b.  DE CLASSIFICATION/DOWNGRADING  SCHEDULE 


4.  PERFORMING  ORGANIZATION  REPORT  NUMBER  (8) 


6*.  NAME  OF  PERFORMING  ORGANIZATION 
Titan  Systems,  Inc. 


Go  ADDRESS  (Cty.  Slut  tnS HP Codt) 

La  Jolla  Gateway,  Suite  500 
9191  Townc  Centre  Drive 
San  Diego,  CA  92122-1231 


6a.  NAME  OF  FUNDINO/8PONSORING  ORGANIZATION 
Defense  Advanced  Research  Project  Agency 


Approved  for  public  release;  distribution  Is  unlimited. 


5,  MONI  TORINO  ORGANIZATION  REPORT  NUMBER (S) 
NOSC  TD  1226 


7a,  NAME  OF  MONITORING  ORGANIZATION 
Naval  Ocean  Systems  Center 


7b.  ADDRES8  (Cty,  Sinned  ZIP  Cede) 

San  Diego,  California  92152-5000 


9.  PROCUREMENT  INSTRUMENT  IDENTIFICATION  NUMBER 
N66001-B4-D-0009 


6c  ADORE  S  8  ( Cly .  Sitlt  me  ZIP  Cede) 

10,  SOURCE  OF  FUNDING  NUMBERS 

PROGRAM  ELEMENT  NO, 

PROJECT  NO, 

TASK  NO, 

Strategic  Technology  Office 

1400  Wilson  Blvd, 

Arlington,  VA  22209 

62301 E 

CM06 

11  TITLE  (nciuet itcutty Ctitilhctlon) 

Daylight  and  Twilight  Sky  Radiance  and  Terrestrial  Irradlance 


12,  PERSONAL  AUThOR(S) 
C.H.  Wahfmsn 


13a,  TYPE  OF  REPORT 


iam 


16,  PAGE  COUNT 

25 


16,  SUPPLEMENTARY  notation 


17,  COSAT I  CODES 


SUB-GROUP 


16.  8UBJ 


MS  I tommuo  on  mttu  1  nnetuty  md  demiy  by  block  numbtr) 


solar  flux 
scattered  flux 
Rayleigh  scattering 


sky  radlsnce 
terrestrial  Irradlance 


19.  ABSTRACT  [Connut  on  mnu  b  ntenwy  me  yStnUy  by  btoebiiuiribei) 

When  solar  radiation  Is  transmitted  through  the  atmosphere  It  may  be  reflected,  transmitted,  or  absorbod,  The  total  effect  of  the 
atmoaphere  upon  a  beam  of  sunlight  passing  through  II  Is  lo  roduce  Its  intensity  by  amounls  varying  with  latitude,  season,  and  degree 
of  cloudiness.  The  atmosphere  weakens  solar  energy  through  1.  molecular  and  small  particle  scattering,  2,  diffuse  reflection  by  larger 
particles,  and  3.  selective  absorption  as  due  to  ozone  and  water  vapor,  The  scattering  and  reflecting  processes  send  a  part  of  the  solar 
energy  back  Into  space,  but  some  of  It  reaches  the  earth's  surface  as  diffuse  sunlight. 

Titan  Systems,  Inc,,  developed  a  model  for  the  sky  radiance  and  terrestrial  Irradlance  and  a  computer  program  was  written  to 
carry  out  the  calculation, 


20.  DlSTRlBUTlON/AVAlUABILITY  OP  ABSTRACT  21.  ABSTRACT  SECURITY  CLASSIFICATION 

□  UNCLASSIFItO/UNUMITtP  [gj  SAMe  AS  RPT  Q  OTIC  USERS  UNCLASSIFIED 


22a.  NAME  OF  RESPONSIBLE  INDIVIDUAL  22b,  TELEPHONE  {MU*  Aim  Cede) 

i,  Mesnan  (619j  553-3737 


220.  OFFICE  SYMBOL 
Code  843 


DD  FORM  1473,  84  JAN 


S3  APR  EDITION  MAY  BE  USED  UNTIL  EXHAUSTED 
ALL  OTHER  EDITIONS  ARE  OBSOltTC 


UNCLASSIFIED 

SECURITY  OLASBUTCATION  OF  THIS  I'AQg 


Daylight  &  Twilight 
Sky  Radiance  &  Terrestrial  Irradiance 

Page  1 


1.  introduction 

When  solar  radiation  is  transmitted  through  the  atmosphere  it  may  De 
reflected,  transmitted,  or  adsorbed.  The  total  effect  of  the  atmosphere  upon  a 
beam  of  sunlight  passing  through  it  is  to  reduce  its  intensity  by  amounts 
varying  with  latitude,  season,  and  degree  of  cloudiness.  The  atmosphere 
weakens  solar  energy  through  (I)  molecular  and  small  particle  scattering,  (2) 
diffuse  reflection  by  larger  particles,  and  (3)  selective  absorption  as  due  to 
ozone  and  water  vapor  The  scattering  and  reflecting  processes  send  a  part  of 
the  solar  energy  back  Into  space,  but  some  of  it  reaches  the  earth's  surface  as 
diffuse  sunlight. 

Titan  Systems,  Inc.  is  interested  to  develop  a  predictive  model  of 
monochromatic  flux  density  reaching  the  earth's  surface.  The  twilight  time, 
when  direct  sunlight  Is  no  longer  available,  is  considered  to  be  of  paramount 
Importance  for  the  model. 

A  model  was  developed  for  the  sky  radiance  and  terrestrial  I'rsdiance  and  a 
computer  program  was  written  to  carry  out  the  calculation.  A  technical 
discussion  of  the  problem,  is  presented  in  Section  2  end  the  computer  program 
is  described  in  Section  3.  The  results  for  a  specified  set  of  conditions  ore 
given  in  Section  4.  Finally,  recommendations  for  improvements  to  the  model 
are  set  forth  in  Section  5.  A  list  of  nomenclature  and  references  are  given  In 
Sections  6  and  7,  respectively. 

2.  Technical  Discussion 
According  to  McCartney  (1975) ... 

'Rayleigh  scattering  by  molecules  represents  the  irreducible  minimum  of 
scattering  along  an  atmospheric  path,  At  lower  altitudes,  file  scattering  by 
particles  always  predominates,  hut  on  average  decreases  more  rapidly  with 
altitude  than  does  Rayleigh  scattering.  This  is  to  be  expected,  because 
representative  scale  heights  of  the  haze  aerosol  are  near  1  km,  while  scale 
heights  of  the  permanent  gas  envelope  are  In  the  range  of  6  to  9  km  and  greater. 
Thus  the  haze  atmosphere  has  a  limited  vertical  extent,  except  for  several 
tenuous  layers  of  particles  such  as  the  one  near  20  km,  In  contrast,  the 
Rayleigh  or  molecular  atmosphere  continues  upwards  to  great  altitudes,  but  Its 
transparency  is  sufficient  that  lines  of  sight  have  lengths  limited  by 
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geometric,  not  optical,  factors.  The  Rayleigh  atmosphere  is  a  useful 
starting  point  for  studying  the  optical  complexities  of  the  real 
atmosphere." 

McCartney’s  analysis  of  the  intensity  of  scattered  flux  will  be  followed 
Although  the  spherical  geometry  of  Figure  1  is  generally  applicable,  for 
simplicity  here  the  earth  and  its  atmospheric  envelope  are  assumed  to  be 
plane-parallel,  as  shown  in  Figure  2.  This  assumption  does  not  cause  errors 
greater  than  2 %  for  all  solar  zenith  angles  less  than  80  deq.  according  to  the 
author.  Clearly,  this  is  not  valid  for  analysis  of  the  twilight  sky.  But  it  is  easy 
to  understand  physically  and  will  provide  the  framework  for  the  twilight  sky 
analysis. 


ZtnitTi 

z 


Figure  1  Geometry  of  atmospheric  scattering  of  sunlight  (after  McCarthy  { 1 9751) 


Considering  Figure  2,  the  value  of  the  sky  radiance  in  direction  OP  is  to  be 
determined.  The  point  P  at  altitude  Z,  lies  within  the  atmosphere,  and  the 
observer  is  at  point  0,  whose  altitude  may  be  taken  as  zero.  (This  is  not 
necessarily  sea  level,  it  may  be  any  reference  altitude)  The  analysis  will  be 
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done  in  terms  of  atmospheric  reduced  height.*  The  reduced  height  of  the 
atmosphere  above  P  is  denoted  Z ,  and  the  total  reduced  height  above  0  by  Z.v 
All  the  solar  rays  incident  on  the  atmosphere  are  considered  parallel,  and  the 
relative  bearing  (or  azimuth)  of  the  line  OP  from  the  line  OS  is  denoted  by  6. 
The  zemth  ang'es  of  the  sun  8nd  the  point  P  are  and  £p ,  respectively  Only  a 
narrow  spectral  interval  of  radiant  flux  is  considered,  and  the  wavelength  is 
taken  to  be  in  a  spectral  region  free  of  absorption. 


Figure  2.  Geometry  of  etmospheric  scattering  of  sunlight  (after  McCartney  [  1 975]). 


♦In  effect,  the  atmosphere  is  replaced  by  a  uniform  density  atmosphere  which 
has  the  same  pressure  and,  in  fact,  total  number  of  molecules  along  a  vertical 
path,  as  the  true  atmosphere,  thus  the  reduced  height  atmosphere  is  finite  in 
extent  The  molecular  optical  thickness  of  a  vertical  path  from  a  given  altitude 
to  the  top  of  the  atmosphere  is  the  product  of  the  reduced  height  for  that 
altitude  and  the  volume  total  scattering  coefficient  relevant  to  the  density  of 
the  reduced  atmosphere  The  reduced  height  at  a  given  altitude  can  be  shown  to 
be  2  =  P/p  g. 
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The  solar  flux,  having  the  irradiance  E?<q  at  the  top  of  the  atmosphere,  is 

attenuated  by  scattering  as  it  travels  the  slant  path  through  the  reduced  height 
Z  to  reach  point  P.  The  irradiance  of  a  unit  volume  at  P  is 

Ek=  EKq exp(-{3m  Z  secCs )  W 

The  scattered  intensity  in  the  direction  P 0  from  a  differential  element  of  path 
at  P  having  a  unit  cross  section  and  a  length 

ds  =  sec  £p  dZ  (21) 

and  characterized  by  the  angular  scattering  coefficient  |3ml0s)  is  given  by 

dlx  =  E^  expl-pm  Z  sec£s )  Pm(e$)  sec£p  dZ'  (3) 

The  attenuation  along  the  path  PO  is  given  by 

exp(-pm(Z0  -  Z  )sec£p)  (4) 

so  that  the  differential  intensity  perceived  at  0  is 

dl\  =  E^  exp(~(3m  Z  sec £s )  pm(0$)  (Z0  -  Z  )  sec £p ]  sec£pdZ  (51 

The  intensity  at  0  due  to  all  elements  8l ong  PO  is  given  by 
Zo 

\K  =  ExoPn(0s)Jexp(-PmZ'sec^)expl-(3m(Zo"Z')sectpl  secCp  dZ'  (6) 
o 

where  it  has  been  implicitly  assumed  that  the  (total  and  angular)  scattering 
coefficients  are  constant  in  the  atmosphere. 

Since  the  total  scattering  coefficient  and  the  zenith  angles  are  independent  of 
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the  altitude,  the  above  integration  can  be  readily  carried  out  in  closed  form. 
First  rewrite  Eq.  (6)  as  follows 


h-  EAofSmteslexpl-proZo  sec£p)sec$p  expl-|3m(sec£s-sec£p)Z  ]  dZ  (7) 


I 


The  integration  yields 


Pml©5* 


Pm 


exp(-pmZ0  secCs)  -  exp(-pm  Z0  sectp, 


(8) 


1  -  sec  £$  /  sec 


for  the  intensity  of  scattered  solar  flux  at  the  bottom  of  a  plane-parallel 
atmosphere  in  a  specific  direction.  Because  the  area  (i.e.,  cross  section)  of  the 
path  was  incorporated  in  the  analysis  this  is  not,  strictly  speaking,  8n 
intensity  (i.e.,  UJ/sr)  but  rather  is  a  radiance  (i.e.,  UJ/m2-sr). 

Equation  (8)  was  derived  for  a  plane-parallel  atmosphere  and  is  not  valid  in 
two  limiting  cases:  observation  along  a  horizontal  line  of  sight  (£p  =  rf/2)  and 
sunset  (ts  =  tt/2).  In  the  former  case  the  sky  radiance  becomes  zero  because 
the  optical  path  becomes  infinite  while  in  the  letter  case  the  sky  radiance 
becomes  zero  immediately  upon  sunset  because  there  is  no  light  at  the  top  of 
the  atmosphere.  Naturally,  these  would  not  be  true  for  a  spherical  atmosphere. 

We  can  formally  derive  an  expression  for  the  sky  radiance  for  a  spherical  earth 
by  accounting  properly  for  the  optical  airmsss  and  solar  irradiance  at  the  top 
of  the  atmosphere.  First  consider  the  optical  air  mass,  m.  For  a  plane-parallel 

atmosphere 


m(£p]  =  secCp  (9) 

while  for  spherical  earth  (refer  to  Figure  3)  we  C3n  use  an  empirical  expression 
given  by  Rozenberg  1 1 9661,  viz., 
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m(£p)  =  (cosCp4  0.025  e'n *  (10) 

Since  this  expression  depends  only  upon  £p  and  is  independent  of  2  we  C8n 
formally  replace  sec£p  in  Eq.  (6)  by 


Next  consider  the  solar  irradiance  of  a  spherical  earth.  McCartney  (19751  has 
considered  the  problem  of  evaluating  terrestrial  paths  at  high  altitudes  and 
zenith  angles  greater  than  about  75°  The  solution  is  provided  by  the  Chapman 
function  which  is  defined  as  follows 
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where  x  =  tZ  ♦  R,J/Hp  and  =  RT/(g  -  Ry)  The  parameter  x  is  seen  to 
depend  upon  the  altitude  of  the  observed  point.  Using  the  Chapman  function  m 
the  derivation  for  the  sky  radiance  would  commence  with  the  solar  irr8diance 
at  P  given  by 


E*  =  exp(-pm  Z  Ch(x,£s) )  (12) 

This  expression  accounts  for  solar  irradiance  on  the  atmosphere  of  a  spherical 
earth  but  does  not  account  for  refraction  of  light  by  the  atmosphere. 

With  Eq.  ( 1 2)  we  can  formally  replace  sec  £s  in  the  analysis  by  Ch(x,£s).  Eq.  (6) 
then  becomes 


^0 

I*  =  EKOpjes)  Jexp(-(3mZ'Cn(x,Csl)  expl-prn(20-Z')m(Cp)]m(Cp)  dZ’  03) 


Tabulated  values  of  the  Chapman  function  show  that  it  is  only  weakly 
dependent  upon  altitude  and  may  be  treated  as  a  constant  as  far  the  integration 
is  concerned.  Then  the  resultant  sky  radiance  is  Identical  to  the  previous 
integration  with  sect,  and  sec£P  replaced  by  Ch(x,C5)  and  m(Cp),  respectively. 


W0S)  f  exp(-f3mZ„'Ch(x,£sI)  -  exp  l-|3m Z0m[£p))  ] 

■x  =  Exo -  !  -  |  (M) 

fs„  [  1  -  Chlx.t.l  /  sec  c„  ] 

A  more  detailed  treatment  of  the  sky  radiance,  including  the  effects  of 

refraction,  is  given  by  Rozenberg  (19661.  Rozenberg's  results  egree  with 
those  obtained  above  when  refraction  is  neglected  in  his  analysts 
Unfortunately,  the  refraction  terms  are  relatively  intractable  and  could  not  be 
readily  adopted  for  a  practical  calculation. 


With  the  sky  radiance  specified  we  con  now  develop  the  monochromatic  flux 
density  or  monochromatic  irradiance  of  radiant  energy.  This  is  defined  as  the 
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normal  component  of  1^  (the  monochromatic  intensity  or  radiance,  1  e  ,  energy 
per  unit  area  per  unit  time  per  unit  solid  angle  per  unit  frequency  interval) 
integrated  over  the  entire  hemispherical  solid  angle  as  shown  in  Figure  4  Liou 
(1580)  gives  the  following  expression  for  the  flux  density 


J I*  cosie)  dcj 

Cm) 


2TT  TT/2 


I^9,<J»l  cos(0)  sin(©)  d©d<t> 


(15) 


where  [9,4>)  are  the  line-of-sig  ,t  zenith  8nd  azimuth  relative  to  the  solar 
vertical,  respectively,  end  are  equivalent  to  (£p,<$)  of  Figure  1 


Figure  4  Illustration  of  solid  angle  and  representation  in  polai  coordinates  (after  liou  ( 1 980)). 


The  problem  will  be  completely  specified  once  we  h8ve  defined  the  angular  and  ' 
total  scattering  coefficients.  McCartney  [1975)  derives  the  angular 
scattering  coefficient  (molecular)  for  the  case  of  Rayleigh  scattering  in  the 
atmosphere,  for  unpolarized  light 

(3m(8s)  =  iz±)l  ( |  4  C0s295  j 
2  N  A4 


(16) 
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where  the  scattering  angle,  8S  can  be  derived  from  the  spherical  geometry 
shown  in  Figure  1  and  is  given  by 

cose,  =  cos  cesCp +  sin£s  sin£p  cos<5  (IT) 

j 

The  total  scattering  coefficient  is  found  by  integrating  the  angular  coefficient 
over  the  spherics!  solid  angle,  thus 

|3m  =  J  (3^9,)  doo  =  8tt3(  n2  -  i  )2  (18) 

J  3  N  A4 

The  scattering  coefficients  in  Eqs.  (17)  and  (18)  exhibit  the  well-known  A4 
behavior  associated  with  Rayleigh  scattering. 
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3  Cojn.nu ter  Program 

A  computer  program  was  written  in  UCSD  Pascal  to  carry  out  the  douDle 
integration  in  Eg  (15)  of  the  monochromatic  radiance  given  in  Eq.  (14)  The 
solution  Fx  =  Fk(Cs)  Is  a  one-parameter  family  depending  only  on  the  solar 
zenith  The  equation  for  the  Chapman  function,  Eq.  (11)  is  integrated  once  for 
each  value  of  the  solar  zenith  (except  that  sec£s  is  used  for  <  75°) 

All  the  integrations  were  earned  out  by  Simpson's  integration  method  This 
method  requires  the  interval  to  be  divided  into  an  even  number  of  steps  The 
general  formula  for  the  integral  is  given  by 


r,  -  1 

n-  2 

2Sfi  ] 

|  Ax/3 

(19) 

j  odd 

j  rvtrt 

The  algorithm  first  chooses  the  number  of  steps,  n  to  be  equal  to  2  and 
automatically  doubles  the  number  of  steps  until  the  result  is  within  a  specified 
tolerance  (1-  10" 5  in  the  present  calculations).  At  each  step  of  the  successive 
approximations  the  function  is  evaluated  only  at  the  odd  positions.  The  sum  for 
the  even  positions  is  obtained  from  the  sum  of  all  previous  interior  positions, 
both  even  and  odd.  The  integration  routine  is  celled  recursively  to  carry  out  the 
double  integration. 

UCSD  Pascal  does  not  permit  a  function  or  procedure  to  be  passed  as  a 
parameter  to  another  routine.  This  problem  was  solved  by  defining  a  TYPE  to 
specify  the  integrand  (eg.  Chapman,  azimuth,  or  zenith  integration). 
Calculation  of  a  single  azimuth  point  requires  the  complete  integration  over 
the  (1  ine-of -sight)  zenith  range  (0,2tt),  each  point  of  which  is  defined  in  terms 
of  the  sky  radiance. 

The  program  is  not  interactive;  operation  of  the  program  requires  no  input  from 
the  user.  The  range  of  solar  zenith  angles  is  specified  in  the  program  and  all 
physical  constants  and  parameters  are  declared  at  the  head  of  the  program. 
Naturally,  the  program  can  be  changed  and  recompiled. 


A  complete,  annotated  program  listing  is  appended  to  this  report. 
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4.  Results 

One  set  of  calculations  of  monochromatic  flux  density  (or  irradiance)  vs  solar 
zenith  was  carried  out  The  solar  zenith  was  varied  from  0°  to  nightfall,  which 
is  defined  in  terms  of  atmospheric  opacity  as  follows 

opacity  =  10  log  (irrediance  at  noon/irradiance)  =  80  decibels 

The  wavelength  and  solar  irradiance  at  the  top  of  the  atmosphere  were  taken 
as  A  =  459.3  nm  and  =  2.00  LU/rrr^nm.  All  physical  parameters  used  in  the 

calculation  were  evaluated  at  sea  level;  the  complete  list  of  parameters  is 
presented  below. 


Quantity 

A 

E*o 

Z. 

n 

N 

X 


value 
459.3  nm 
2  00  Ul/m2*nm 

8.44  km 
1 000293 
2.68719  •  !025  m*3 
612.5 


The  results  showed  that  under  these  conditions  night  falls  at  a  solar  zenith  of 
ts  =  103”,  or  about  52  minutes  after  sunset.  This  result  is  predicated  on  the 
assumption  that  the  solar  zenith  is  zero  at  the  solar  noon. 

The  complete  tabulated  results  for  Chapman  function,  flux  density,  and  opacity 
are  appended  to  this  report 
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5.  Recommendations 

There  are  several  improvements  which  could  be  made  to  the  model  described  in 
Section  2.  Some  of  these  are  noted  below: 

Phenomenon _  References _ 


refraction  of  light 

McCartney  [19751 
Ro2enberg  [1966! 

shadow  of  earth  on  atmosphere 

Rozenberg  [19661 

airglow 

Henderson  [1977] 
Fleagle  [19631 

selective  absorption 

Liou  [1980] 

particle  (haze)  scattering 

McCartney  [19751 

V  '  - 

multiple  scattering 

McCartney  [19751 

diffuse  scattering  by  clouds 

Liou  [1980] 
McCartney  [19751 

Additionally,  the  program  can  be  adapted  to  permit  calculation  of  the  solar 
2enith  for  a  given  latitude,  day  of  the  year,  and  time  of  day  which  would  permit 
a  more  realistic  calculation  for  a  particular  locale.  Figure  5  shows  the 
geometry  for  calculating  the  solar  zenith.  Ltou  (i960)  gives  this  result 

cos  0S  =  sin  A  sin  <5  ♦  cos  A  cos  5  cos  ft,  (20) 


where 


9S  =  solar  zenith 

7\  =  latitude  of  the  observation  point 

ft,  =  hour  angle  (this  is  zero  at  solar  noon  and  increases 
15°  for  every  hour  before  or  after  solar  noon) 

=  solar  inclination 
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The  solar  inclination  is  a  function  of  the  day  of  the  year  only,  and  is 
independent  of  the  location  of  the  observation  point  It  varies  from  23*27’  on 
June  21  (summer  solstice)  to  -23*27’  on  December  22  (winter  solstice). 
Kreider  and  Kreith  !  1981]  give  the  following  empirical  equation  for  the  solar 
inclination 


sin <5  =  0.39795  cos(0  985631 N  -  173)) 


(21) 


where  N  =  numoer  of  the  day  (January  1st  =  1). 
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Nomenclature 


Sumbol 

Definition 

Units 

E* 

monochromatic  irradiance 

UU/m2-nm 

E*o 

at  the  top  of  the  atmosphere 

tU/m2-nm 

F* 

monochromatic  flux  denity  (irradiance) 

UU/m2 

g 

acceleration  of  gravity 

m2/sec2 

He 

density  scale  height  of  molec.  atmos. 

km 

h 

monochromatic  radiance 

lll/m2-nm-sr 

m 

optical  air  mass 

n 

refractive  index 

N 

number  of  motecules  per  unit  volume 

m'3 

R 

gas  constant 

m2/sec2-°K 

R. 

mean  radius  of  the  earth 

km 

s 

distance 

km 

T 

temperature 

*K 

2 

altitude 

km 

Z 

reduced  height  of  atmos.  (ref.  altitude) 

Zo 

reduced  height  of  atmos.  (sea  level) 

Pmle*) 

volume  angular  scattering  coefficient 
(molecular),  unpolarized  light 

km'1  sr'1 

Pm 

volume  total  scattering  coefficient 
(molecular) 

km'1 

atmos.  lapse  rate  of  temperature 

°K  km'1 

azimuth  relative  to  solar  vertical 

rad 

zenith 

rad 

zenith  of  line  of  sight  at  observer 

rad 

zenith  of  sun  at  observer 

rad 

es 

observation  angle  of  scattering 

rad 

A 

wavelength 

nm 

X 

parameter  of  the  Chapman  function 

GO 

solid  angle 

sr 
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'  K Y  R A  D I  ANT  F  L.  U> 


Monday,  14  January  1935  4:09  PM 


( **«*«*****«*******«***********#«***t««»»**4.»4*, ***,««**, 
* 

*  PR0GR4M  TITLE:  Calculation  of  Sky  Radiant  Flu* 

*  URITTEN  BY:  Cye  H.  Uaicman,  Ph.D. 


CATE  URITTEN: 
URITTEN  FOR: 


8  January  1985 


Gary  Lee ,  Ph.D., 


PROGRAM  ELMttARY:  This  program  calculates  the  monochromatic 
■flux  density  (or  irnadiance)  of  radiant 
energy  by  integrating  the  monochromatic 
intensity  on  radiance"  I  (power/area«+rea*sr) 
over  tne  snnerical  solid  angle  encompassing 
the  sky  hemi sonere . 


1NPIT  &  OUTPUT  FILES: 
SKY. FILE. TEXT 


description 

Output  File  of  tabulated  results 


(  Notes  on  the  integration  scheme:  A  double  Simpson's  integration  is 
perFormpo  over  the  azimuth  and  zenith  angles.  The  integration 
procedure  is  called  recursively  in  order  to  do  this.  Since  UCSD 
Pascai  doesn't  permit  a  FUNCTION  or  PROCEDURE  to  be  passed  as  a 
parameter  to  another  routine  we  circumvent  this  problem  by  using 
a  detmed  TYPE  to  specify  tne  integrand.  When  the  Simpson's  procedure 
calls  the  Function  For  the  integrand  it  chooses  tne  correct  Function 
in  a  CASE  statement.  Calculation  ot  an  azimuth  point  recumes  the 
complete  integration  over  tne  (line  ot  sight)  zenith,  hence  that 
'Function*  is,"  in  Fact,  another  S  inpson' s  ’  i  n  tegrat  i  on  oF  the  intensity.  ) 

(xL  T ! TAN : SKY. LST. TEXT) 

PROGRAM  SkyRadi antFl ux ; 

CONST  {  solar  and  terrestrial  properties  ) 


1 ambda 
E  zero 


459.3; 

2.00; 


2  zero  =  8.44; 
n~=  1.000293; 

N  mol  =  2.6871 9E25; 
cFTi  =  612.5; 
night_deF  *  S0.0; 


{  wavelength  (blue/green) ,  nm  ) 

{  ex tr aterrestr i al  solar  irradiante  at 
Frequency  lambda,  U/(m’2*nm)  ) 

(  reduced  height  oF  atmosphere/sea  level,  km  ) 
<  retractive  index  oF  air  (STP)  ) 

{  loschmidt's  number  (STP),  mol ecul es/m‘3  ) 

{  parameter  oF  the  Chapman  Function  ) 

(  def  of  night  opacity,  decibels  ) 


{  miscellaneous  and  system  constants  ) 


pi  =  3.141593; 
convert  =  1.745329E-2; 
tol  =  1.0E-5; 
bell  =  7; 

FF  =  12; 


debugg 
outF  i  T 


ing  =  FALSE:  (  turns  de 

ename  =  ' TITAN: SKY. F] LE, TEXT 


(  a  transcendental  number  ) 

(  decrees  to  radians  ) 

{  tolerance  For  numerical  integration  } 
(  ASCII  character  for  bell  sound  ) 

{  ASCII  character  For  Form  Feed  } 

(  turns  debugging  on/off  } 

/cue  TCVT'  . 


{  1 imi ts  of  integration  ) 
az  1st  *  0.0:  { 
azHast  =  6.283185;  C 
zen  1st  =  0.0;  ( 
zen~l ast  =  1 .570796;  ( 


{  1 ower  limit  azimuth  ) 

{  upper  limit  azimuth,  2*pi 


{  lower  limit  zenith  (line  of  sicht)  } 

{  upper  limit  zenith  (line  of  sight),  pi/2  ) 


{  descriptions  of  integrands  For  Simpson's  integration  ) 
integrand  =  (chapman,  azimuth,  zenith); 

{  mout/output  parameters,  input  is  1  line  in  SKY. DATA. TEXT  ) 
zeta  S,  C  zenith  of  sun  at  observer,  radians; 

this  could  also  be  calculated  From  the 
latitude,  day  of  year,  &  time  of  day  ) 
Flux:  REAL;  <  result:  monochromatic  Flux  density  ) 


{  independent  variables  ) 


zeta_P, 
del ta:  REAL; 


(  zenith  of  line  of  sight  at  observer,  rad  ) 
(  azimuth  relative  to  solar  vertical,  rad  ) 


(  scattering  and  refraction  variables  ) 


theta  s, 
beta, 

beta_total , 
chaoman_f  n , 
be  ta_Z : "REAL ; 


(  scattering  angle,  radians  } 

I  angular  scattering  ccet* ,  l/;km*sr)  ) 

(  total  scattering  coefficient.  1/km  ) 

(  Chapman  function  c a  1  c ' d  for  chi  !<  :eta  $  ) 
{  calc'd  product  of  oeta  total  *  Z  zero  T 


(  miscellaneous  variables  ) 


opac i ty , 


(  relative  irradiance  (noon  ref),  decibels; 
this  is  the  relative  capacity  of  matter  to 
obstruct  the  trans-m  i  ss  i  on  cf  radiant  eneroy 
(  monochrcmat  i  c  raoiance  at  noon,  ij/m‘2*rm  .1 


noon_flux:  REAL; 
cnoice:  inteqrand; 
i:  INTEGER;  ‘ 

Chap^quad.exp  S:  DEAL; 
infiTe,  outfiTe:  TE>'.T ; 
separator! me, spaces:  STRING; 

PROCEDURE  Integrate(choice:  integrand;  ! ower .upper , tol :  REAL: 

CAR  sum: REAL); FORWARD; 

FUNCTION  EXP0( x: REAL) : REAL; 

(  checks  for  exponent  to  be  within  bounds/max  real  =  10*308  ) 

BEGIN 

IF  x  <  -700  THEN 
EXPO  :=  0.0 
ELSE 

£<?0  :=  EXP(x); 
e^'D  (EXPO); 

FUNCTION  Scatter i ngFn( sol ar_zen i th ,  obs_pt_:en i th ,  azimuth:  REAL): REAL; 

{  calculates  the  ratio  of  angular  to  total "scat ter i ng  function  as  a 
function  of  the  scattering  ancle  ) 

CAR  Deta_rat i o,cos_of_theta :  REAL; 

BEGIN 

cos_of  theta  :=  COStsolar  zenith)  *  COS(obs_pt_zem th) 

♦  SIN( so! ar“zen i th)  ♦  SIN(obs  pt  zenith)  *  C0S(a: imuth) ; 
beta  ratio  :=  (1  ♦  SQR(cos  of  theta))  »  3  /  (8"*  pi); 

IF  debuoQi nQ  THEN  BEGIN  “  " 

WRITELN( 'C0S( theta  s)  ',cos  of  theta): 

URITELN( 'beta  ratio  =  '.beta  ratio); 

END: 

ScatteringFn  :=  beta  ratio; 

END  {ScatteringFn); 


FUNCTION  ChaoFXtchi .zeta.psi :  REAL). ’REAL: 

(  calculates  the  kernel  tor  the  Chapman  function  integration  ) 
VAR  kernel :REAL; 


BEGIN 

IF  psi  =  0  THEN 
kernel  :=  0 
ELSE 

kernel  :=  EXPOtchi  »  (1 
ChapFX  :*  kernel ; 

END  (ChapFX); 


-  SJN(zeta)  /  SJN(psi)))  /  SQR(SIN(psi )) ; 


FUNCTION  Intensi ty(del ta,  zeta_P:  REAL):REAL; 

{  calculates  the  intensity  function  already  inteorated  in  closed  for  a 
plane  parallel  atmosphere  with  curvature  corrections  for  air  mass  and 
solar  irradiance  ) 

VAR  kernel,  radiance,  airjnass,  exp _P,  cos_zeta:  REAL; 

BEGIN 

cos_:eta  :=  C0S(:eta_P); 

air_mass  ;=  I  /  (cos  zeta  ♦  0.025  *  EXP(-il  *  cos  zeta)); 
exp _P  :=  EXPO(-beta_i  *  air  mass); 

IF  air_mass  =  chapman_fn  THEN 
radiance  :=  E_zero  *  Scat  ter i naFn( zeta_S,  zeta  P,  delta) 

*  beta_Z  »  chapman  fn  *  exp  S 

ELSE 

radiance  :=  E_:ero  *  Scatter i naFn(zeta_S,  zeta_P,  delta) 

*  (exp  s  -  exp_PJ  /  (1  -“chapman  fn/airjnass) ; 

(  multiply  by  SIN  *  COs  to  oet“kerne>  ,'or  flux  density  Integration  ) 
kernel  :*  radiance  *  C0S(:efa  P)  *  SIN(zeta  P) ; 

Intensi ty  :*  kernel ; 

END  (Intensity); 
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FLECTION  FXtchcice:  integrand;  x:  REAL):REAL; 

(  "choice1  gives  the  correct  function  for  Simpson's  integration  ) 
VAR  result:  REAL; 

BEGIN 

CASE  choice  OF 
azimuth:  BEGIN 

cel ta  :=  x ; 

Integratevzem th,  :en  1st,  :en  last,  tol  .  result); 
END;  ' 

zenith:  BEGIN 

:eta_?  :=  x; 

IF  debugging  THEN  BEGIN 
LiR!T£Ln(  oe  I  ta  =  '.delta.1: 

UR!TELN('zeta  P  =  '  ,:eta  ?■  ; 

END; 

result  :=  Intensi ty(del  ta,  zeta  P>; 

END; 

chapman:  result  :=  ChapFX(chi,  zeta  S,  x); 

END  (of  CASEs); 

FX  :=  resu 1 1 ; 

END  (FX); 


PROCEDURE  Inteqratetchoice:  integrand;  1 ewer  .upper , tol :  REAL:  VAR  sum:  REAL); 
(  numerical  integration  by  Simpson's  method  of  FX(x)  for  x  =  (lower, upper]  } 
VAR  , .steps: INTEGER; 

x , de 1 1  a_x , euen_sum , ood_sum , end_sum , suml : REAL ; 

BEGIN 

x  :=  lower; 
steps  :=  2; 

delta_x  :=  (upper  -  lower)  /  steps; 
cdd_sum  :=■  -/(choice,  lower  ♦  delt«_x>; 
euen_sum  :=  0; 

end_sum  :=  FX(choice,  lower)  ♦  RXtchoice,  upper;; 
sum_:=  tend  sum  ♦  4  *  odd_sum>  »  delta_x  /  3; 

IF  debugging  THEN 
BEGIN 

CASE  choice  OF 

azimuth:  URlT£LN(step$:6,  x : 1 0 : 3 ,  spaces,  sum:iC:3); 
zenith:  URITELN(spaces,  steps:©,  x  :  1 0 : 3 ,  sum:i0:3>; 

END  (of  CASEs); 

END; 

REPEAT 

steps  :=  steps  *  2; 
suml  :=  sum; 

delta_x  :=  (upper  -  lower;  /  steps; 
even_sum  ewen_sum  ♦  odd_sum; 
odd  sum  :=  0: 

FOR-.  :=  1  TO  steps  DIV  2  DO 
BEGIN 

x  ;=  lower  ♦  delta_x  *  (2  *  i  -  1); 
odd  sum  :=  odd  sum  *  FX( choice,  x); 

END;  ~ 

sum  :=  (end_sum  ♦  4  »  odd_sum  ♦  2  *  ever>_sum)  *  delta_x  /  3; 

IF  debugging  THEN 
BEGIN 

CASE  choice  OF 

azimuth:  WRITELN<steps:6,  x : 1 0 : 3 ,  spaces,  sum:l0:3); 
zenith:  URITELNt spaces,  steps:6,  x : 1 0 : 3 ,  sum : 1 0 : 3 ) ; 

END  (of  CASES); 

END; 

UNTIL  A6S(sum  -  suml)  <=  A6S(to!  *  sum); 

END  (Integrate); 


PROCEDURE  Initialize; 

(  various  and  sundry  business  to  take  care  of  ) 


BEGIN 

(  calc  total  scattering  coefficient  in  1/km  ) 
beta  total  :=  9  *  pi  *'S0R(pi)  *  SQR(SQR(n;-l )  *  1 . 0E39  / 
<3~*  N  mol  *  SQR< SQR( 1 ambda) ) ) ; 

(  calc  tfie  Chapman  function  for  the  given  solar  zenith  ) 
IF  zeta_S  <  1.3  (74.5  degrees)  THEN 
chapman  fn  :=  I  /  C0S(zeta  S) 

ELSE 

BEGIN 

lnteorate(chapman,  0,  zeta  S,  tol,  Chap_5uad>; 
chapman^fn  :=  chi  *  SlN(zefa_S)  *  Chap_quad; 
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END: 

beta  l  :=_beta^total  *  2  zero: 
exp  5  :=  t\PQ<  -be  t  •-  Z  *  chapman  -f n ) ; 

IF  CebugQing  THEN  _E6IN 
WRITECNt ' be ta_to cal  =  ' ,beta_total ) ; 
URITELNC ' chapman  *n  =  '.chapman  tn); 
URITELNC 'beta  2  =  ',beta  2); 
WRITELN('exp  p  =  '  ,exo  3); 

END: 

separator t i ne  :=  ; 

spaces  :=  '  ' ; 

FOR  i  :=  1  TO  7?  DO 

separator ]  i  ne  :=  CCNCATC separator 1 i ne ; 
END  (Initial i ze } ; 


PROCEDURE  PaqeHeader; 

(  write  input  data  to  output  -file  and  titles  tor  tabulated  output  > 

BEGIN 

WRITELNdsutf  Me):  URITELNC  ou  t  f .  1  e  ) ;  URITELNC  ou  t  f  i  1  e  ) ; 

WRITELNC ou tf i 1 e , ;Ca! cu 1  at i on  of  Sky  Radiant  Flux'1*; 

URITELNC ou  1 f  i  1  e) :  URITELN(outf  Me); 

WRITE  <  ou  t  f  1 1  e  , '  Sol  ar  Zem  th' :  15,  'Chapman  ':18,'Mono.  Flux  Density  :  22 )  ; 
URITELNCoutf  Me, '  Opacity  ' :  1 5) : 

WRITE  <  ou  t  -f  1 1  e  , '  (decrees)  ' : 15, 'Func t i on' : 18, '  (Wat ts/m‘ 2*nm)  ':22); 

WR I TELN( ou 1 1 1 1 e . ' ( dec i be  1 s ) ' : 1 5) ; 

WRITE  ( ou  1 1 1 1  e  , ' - '  :  1 5. ' - ';!8,' -  :  22) : 

WRITELN( ou tf  Me,' . . '  :15> ; 

(  and  write  it  cut  the  the  screen  as  well...  ) 

WRITELN:  WRITELN;  WRITELN: 

URITELNC ' Cal cul at i on  of  Sky  Radiant  Flux'); 

WRITELN:  WRITELN: 

WRITE  ('Solar  Zen i th' : 15, 'Chapman  ':18,'Mono.  Flux  Densi ty' :22) ; 
WR!TELN('  Ooacitv  'M5-; 

WRITE  <’  (decrees;  ' : 15, 'runct i on' : 18, '  (Uat ts/m‘ 2*nm)  ':22); 

WRITELN( ' (decibels)' : 1 5) : 

WRITE  (' - ' :  1 5, ' . 'M3,' . — '  :  22) ; 

WRITELN( ' - ' :  1 5) ; 

END  (PageHeader ) j 


BEGIN  (SkyRadi antFlux) 

REUR1TEC  ou  tf  Me,  out*  i  1  ename) ; 

PaaeHeader : 
zeta_S  :=  0: 
opacity  :=  0; 

WHILE  opacity  <=  night_def  DO  BEGIN 
Initialize; 

IntegrateCaz imuth ,  az_lst.  az_last,  tol ,  flux); 

IF  zeta  S  =  0  THEN  noon_flux  :=  flux; 
opacity-:*  10  *  L06Cnoon  flux/flux); 

WRlTELNCoutf  M  e ,ze ta  S/convert :  15:2, chapmr.n  f  n :  18: 3, f  1  ux  : 22 : 9 , 

opac i ty: 15:2) ; 

URITELNC zeta  S/con vert : 15:2,chaoman  f n : 1 8 : 3 , f lux :22:9,opaci ty: 15:2) ; 
zeta  S  :=  zeTa  S  ♦  1  *  convert; 

END  (WHILE): 

CLOSECoutf  Me,  LOCK); 

END  (SKyRadi antFl ux) . 


M ON OCHR DMA T I 


Z  EM  I  TH 
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E  (Solar  irradiance;  =  2.00  u/m‘ 

'2*nm  at  (wavelength) 

=  459.3  nm 

Sol ar  Zenith 

Chaoman 

Mono .  Flux  Dens- 1  ty 

Ooac 1 1  v 

degrees-* 

Function 

(U!atts/m‘2*nm> 

t dec i ce 1 s ) 

0.00 

1 .000 

0.277675551 

0.00 

5.00 

1  .004 

0.277429912 

0.00 

10.00 

1.015 

0 . 276693685 

0.02 

15.00 

1  .035 

0.275468486 

0.03 

20.00 

1 .064 

0 . 273755427 

0.06 

^5.00 

1.103 

0.271552579 

0.10 

30.00 

1.155 

0.268350962 

0.14 

35.00 

1  .221 

0 . 2 c>5628529 

0.19 

40.00 

1  .305 

0.261841014 

.  0.25 

<15.00 

1  .414 

0.257407733 

0.23 

50.00 

1  .556 

0.252188539 

0.42 

55.00 

1.743 

0.245944344 

0.53 

60.00 

2.000 

0.238264905 

0.66 

65.00 

2.366 

0.22S426001 

0.85 

70.00 

2.924 

0.215079046 

1.11 

75.00 

3.782 

0.197046060 

1  .49 

76.00 

4.033 

0.192224544 

1  .60 

77.00 

4.321 

0.186924108 

1  .72 

78.00 

4.653 

0.131072370 

1  .£6 

Z5, .  00 

5.040 

0  •  i  <  4  f  S  6  4  7  6 

2.02 

80.00 

5.496 

C . 1 67366458 

2.20 

81  .00 

6.040 

0.159299163 

2.41 

82.00 

6.700 

0.150256251 

2.67 

83.00 

7.515 

0.140098668 

2.97 

84.00 

8.539 

0.123688104 

3.34 

85.00 

9.861 

0.115911849 

3.79 

86.00 

11.616 

0.101732926 

4.36 

87.00 

14.027 

0.086277843 

5.03 

88.00 

17.491 

0.069965942 

5.99 

89.00 

22.689 

0.053632170 

7.14 

90.00 

31 .037 

0.038481520 

8.58 

91.00 

45.449 

0.025662122 

10.34 

92.00 

72.638 

0.015735736 

12.47 

93.00 

129.575 

0.008696670 

15.04 

94.00 

264.030 

0.004230336 

18.17 

95.00 

627.384 

0.001771989 

21 .95 

96.00 

1765.305 

0.000628543 

26.45 

97.00 

5935.977 

0.000186841 

31.72 

98.00 

23955.375 

0.000046309 

37.73 

99.00 

116186.717 

0.000009553 

44.63 

100.00 

677317.886 

0.000001640 

52.29 

101  .00 

4743964.314 

0.000000234 

60.74 

102.00 

39898639.363 

0.000000023 

69.98 

103.00 

402682251.982 

0.000000003 

S0.02 

